%%Run with sbatch -p ncf_holy -n 1 -t 12-12:00 --mem=30000 --wrap="matlab -nodisplay -nosplash -nojvm -r $'searchlight_hv(subj#);exit'"
function cf_classification_hadvdid_masks(sub,val)

model='B81'
%classifier = 'gnb_MATLAB';%assumes no correlation among voxels
classifier = 'svm_MATLAB';
%classifier = 'lda_shrinkage';%assumes some correlation among voxels and that that correlation is the same across conditions
%classifier = 'svm_linear'
ROIs=15;

if sub < 10
    s = ['0',num2str(sub)]
else s = num2str(sub)
end

cd (['~/../../ncf/cushman/Regan/desert_island/desisl_',s,'/desisl_',s,'.results']);

%make a directory to move all of the output to
if ~exist('hadvdid_class','dir')
    unix(sprintf('mkdir hadvdid_class'));
end

%%
%1. Load data
for curroi=1:ROIs

    if curroi==ROIs
        r='full_mask'
    else
        r=['roi',num2str(curroi)]
    end

    %Load brain data
    load(['s',s,'_BLOCK81_LSA_braindata_hadvdid_',r,'.mat'])

    %Load folds & regressors
    load (['~/../../ncf/cushman/Regan/desert_island/desisl_',s,'/desisl_',s,'.results/s',s,'_mvpa_params.mat']);

    if val==1
        labels=lv_labels;
        conds=lv_conds;
        obj='lv'
    elseif val==2
        labels=hv_labels;
        conds=hv_conds;
        obj='hv'
    end
      
    %orient labels and conds with data and runs    
    labels=labels';
    conds=conds';

    sizelabels = size(labels)
    sizedata = size(data)

    %create new folds array to reflect the fact that each condition only occurs
    %every other trial
    ctr=1;
    for run_num=1:8
        cut_folds(ctr:ctr+3)=run_num;
        ctr=ctr+4;
    end

%%
%3. Specify which labels, folds, and data will be used for this classificaiton depending on
%the condition
for cur_cond=1:max(conds)



    if cur_cond==1
        c='targ_had'
    else if cur_cond==2
            c='targ_did'
        else if cur_cond==3
                c='third_had'
            else
                c='third_did'
            end
        end
    end

    cut_labels=labels(conds==cur_cond);
    cut_data=data(:,(conds==cur_cond));

    sizecutlabels=size(cut_labels)
    sizecutdata=size(cut_data)
    sizecutfolds=size(cut_folds)

    nFolds=length(unique(cut_folds))
    [nVoxels,nTrials]=size(cut_data)

%%
%4. Run classificaiton

%Create results arrays
    accuracy=[];
    full_scores=[];
    full_guesses=[];

    for cur_fold=1:nFolds
    
        %Output the current fold for updating
        cur_fold
    
        %Create more results arrays
        cur_run_scores=[];
        cur_run_guesses=[];
    
        %Specify which labels are training and which are test labels
        test_labels=cut_labels(cut_folds==cur_fold);
        train_labels=cut_labels(cut_folds~=cur_fold);

        %Identify training and test data
        test_data=cut_data(:,cut_folds==cur_fold);
        train_data=cut_data(:,cut_folds~=cur_fold); 
            
        train_data=train_data';
        test_data=test_data';
            
        %Train classifier with training data
        Mdl=fitcsvm(train_data,train_labels);
         
        %Test classifier on test data
        [test_output,score]=predict(Mdl,test_data);
         
        %Collect the classification score for each trial
        for cur_trial=1:length(test_labels)    
            %Get the classification score for each trial for the correct
            %answer, where positive values means it is right and negative values means it is wrong:
            cur_run_scores(cur_trial)=score(cur_trial,test_labels(cur_trial));       
        end
        
        %Compute the average score across all trials in this run
        curroi_avg_scores(cur_fold)=mean(cur_run_scores);
            
        %Compare test classifications to actual label
        test_output=test_output';
        correct=(test_output==test_labels);
        cur_run_guesses=[cur_run_guesses;test_output];
         
        %Comput the proportion correct creating a matrix with a column for
        %each run
        accuracy(cur_fold)=sum(correct)/length(correct);
 
    
        %full scores and full guesses are a matrix of the scores and
        %guesses across runs, with a row for each run and a column for each
        %trial in that run
        full_scores=[full_scores;cur_run_scores];
        full_guesses=[full_guesses;cur_run_guesses];
    
        %save posts and accuracies through this run
        save(['~/../../ncf/cushman/Regan/desert_island/desisl_',s,'/desisl_',s,'.results/',model,'_',obj,'_',c,'_full_scores_',classifier,'_hadvdid_',r,'.mat'],'full_scores')
        save(['~/../../ncf/cushman/Regan/desert_island/desisl_',s,'/desisl_',s,'.results/',model,'_',obj,'_',c,'_full_guesses_',classifier,'_hadvdid_',r,'.mat'],'full_guesses')
        save(['~/../../ncf/cushman/Regan/desert_island/desisl_',s,'/desisl_',s,'.results/',model,'_',obj,'_',c,'_accuracies_',classifier,'_hadvdid_',r,'.mat'],'accuracy')
end

%average together the accuracy and scores for each ROI across runs, so there is a
%is a row for each roi and a column for each condition
avg_accuracy(curroi,cur_cond)=mean(accuracy);
avg_scores(curroi,cur_cond)=mean(curroi_avg_scores);

end
end
%save accuracies and scores
save(['~/../../ncf/cushman/Regan/desert_island/desisl_',s,'/desisl_',s,'.results/',model,'_',obj,'_',classifier,'_hadvdid_rois.mat'],'avg_accuracy', 'avg_scores')

%move everything to the directory for this classificaiton
unix(sprintf(['mv B81_',obj,'*mask* ./hadvdid_class']))
unix(sprintf(['mv B81_',obj,'*roi* ./hadvdid_class']))

cd ../../
end
